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' Tricritical point in charge-order systems and its criticality are studied for a microscopic 

£SJ , model by using the mean-field approximation and exchange Monte Carlo method in the clas- 

sical limit as well as by using the Hartree-Fock approximation for the quantum model. We 
study the extended Hubbard model and show that the tricritical point emerges as an end- 
point of the first-order transition line between the disordered phase and the charge-ordered 
phase at finite temperatures. Strong divergences of several fluctuations at zero wavenumber 
| are found and analyzed around the tricritical point. Especially, the charge susceptibility Xc 

■ and the susceptibility of the next-nearest-neighbor correlation \R are shown to diverge and 
their critical exponents are derived to be the same as the criticality of the susceptibility of 
the double occupancy Xd(o)- The singularity of conductivity at the tricritical point is clarified. 

■ We show that the singularity of the conductivity a is governed by that of the carrier density 
and is given as \a — a c \ ~ \g — g c \ pt (A log \g — g c \ + B), where g is the effective interaction 
of the Hubbard model, a c (g c ) represents the critical conductivity(interaction) and A and B 
are constants, respectively. Here, in the canonical ensemble, we obtain p t — 2/3 t = f/2 at the 
tricritical point. We also show that pt changes into p' t = 2f) = 1 at the tricritical point in the 
grand-canonical ensemble when the tricritical point in the canonical ensemble is involved within 
the phase separation region. The results are compared with available experimental results of 

>-^j ' organic conductor (DI-DCNQI)2Ag. 

KEYWORDS: tricritical point, extended Hubbard model, charge order, organic conductor, charge suscep- 
O , tibility, doublon susceptibility, phase separation 
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1. Introduction 

Charge orderings and charge density waves are widely 
observed in various correlated electron systems such as 
transition metal compounds 1 and organic conductors. 2 
Such regular alignment of electrons with a periodicity 
longer than that of the unit cell at high temperatures 
is stably realized, particularly when the electron den- 
sity is commensurate, where the number of electrons per 
unit cell is a simple fraction. 3 These ubiquitous phenom- 
ena have attracted recent intensive interest not only be- 
cause of its own right but also because dramatic phenom- 
ena such as colossal magnetoresistance in perovskite-type 
manganese oxides are observed immediately when the 
charge order melts. Charge orderings are in many cases 
competing or coexisting with magnetism, ferroelectric- 
ity and superconductivity, resulting in strong coupling to 
transport and magnetic properties. Charge orderings are 
in fact frequently accompanied by metal-insulator transi- 
tions through the growth of charge-order parameter and 
the opening of gaps at the Fermi surface. 

Sensitive changes in transport, optical, dielectric and 
magnetic properties at or near charge-order transi- 
tions have stimulated intensive research on their con- 
trol through the charge-order transitions. To understand 
various possibilities, it is desired to clarify the basis of 
charge-order transitions themselves, since the sensitive 
changes of physical properties and emergence of compet- 
ing phases may be deeply influenced by the nature of the 
transition. 

Charge-order transitions take place either as first-order 
or continuous ones. They even show both continuous 



and first-order boundaries meeting at a tricritical point 
(TCP) in phase diagrams as in (DI-DCNQI) 2 Ag 4 (see 
Fig. 1). TCP is characterized as an endpoint of the first- 
order transition line. However, the first-order line does 
not simply terminate as the normal critical point because 
the charge-ordered phase explicitly breaks translational 
symmetry and the phase boundary should not disappear. 
After the termination of the first-order line, it continues 
as the continuous transition line (the lambda line) and 
the meeting point of the first-order and the lambda lines 
is called TCP. 5 For (Dl-DCNQl) 2 Ag, the transition to 
the charge-ordered phase actually takes place with a de- 
crease in either temperature or pressure as we see in Fig. 
1. This means that the charge-order transition may be 
caused not only by suppressing thermal fluctuations but 
also as a quantum transition by suppressing the band- 
width with decreasing pressure, ft offers an intriguing 
field of quantum critical phenomena for second-order and 
tricritical transitions. 

In this paper, we focus on various fluctuations ex- 
pected around the continuous boundary (namely the 
lambda line) as well as around the TCP. To the authors' 
knowledge, there exist no systematic studies on the ques- 
tion which fluctuations and their divergences character- 
ize the criticalities of the charge-order transitions partic- 
ularly around the TCP. Divergences of fluctuations and 
singularities of physical properties are the main subjects 
of this paper. We will show that for charge-order tricrit- 
icality, two completely different fluctuations diverge: As 
is well known, one is the order parameter susceptibility 
at the ordering wavenumber, which leads to the Bragg 
scattering of the charge density response in the ordered 
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Fig. 1. Schematic phase diagram of (DI-DCNQI)2 Ag. 

phase. The other divergence occurs in some of density 
fluctuations at vanishing wavenumber. This divergence 
also triggers the singularity of specific heat and further 
nontrivial effects arise. We will also show how the singu- 
larity of conductivity appears at the critical and tricriti- 
cal points. At TCP, we will show that the singularity of 
conductivity in the canonical ensemble is different from 
that in the grand-canonical ensemble. Available experi- 
mental results on the conductivity are compared with the 
results of present study. The experimental results could 
be interpreted by the mean-field criticality of TCP in 
the grand-canonical ensemble. We then propose further 
experiment to reveal the true tricriticality. 

In the Ginzburg-Landau (GL) scheme, TCP is ex- 
pressed by the (j) 6 theory. 5 The free energy F is expanded 
up to the sixth order with respect to the order parameter 
m as 

F = rm + um + vm 6 . (1) 

If u > 0, r = represents a conventional Ising-like critical 
point. If u < 0, three minimum states can emerge. These 
three minima represent the coexistence of three phases. 
Then three phase boundaries appear between each com- 
bination of phases as first-order phase boundaries. Since 
the coexistence of two phases is represented by a surface 
in the parameter space, the coexistence of three phases 
is represented by a single curve, whose endpoint is TCP. 
In other words, TCP appears as the crossing point of 
three lambda lines, where the lambda lines are termi- 
nating lines of each first-order phase-boundary surface. 
Figure 2 shows a schematic phase diagram of TCP in the 
space, where T, g and W represent temperature, 
interaction and the field which is conjugate to the order 
parameter, respectively. Near TCP, r and u are described 
by the linear combination of \T — T c \ and \g — g c \, where 
T c and g c represent the critical temperature of TCP and 
the critical interaction of TCP, respectively. The first- 
order surfaces are depicted as shaded surfaces. In eq. 



Lambda line 




Fig. 2. Schematic phase diagram of TCP in the T-g-H^ space. T, 
g and H" 1 represent temperature, interaction and the field which 
is conjugate to the order parameter, respectively. The shaded 
surfaces represent the first-order transition surfaces. 

(1), u = and r = represent TCP. TCP belongs to 
a different universality class from that of the continuous 
transition lines. It is also known that the upper critical 
dimension of TCP is three so that fluctuations become 
irrelevant for d > 3 and marginal at d = 3. Therefore, the 
mean-field treatments are more or less justified in three 
dimensions except for possible logarithmic corrections in 
the exponents. 

The GL free energy given in eq. (1) describes the tran- 
sition of the scalar order parameter leading to Ising-type 
symmetry breaking. The commensurate charge ordering 
occurring at a simple fractional filling of electrons is in- 
deed represented by a discrete symmetry breaking of spa- 
tial translation. For example, at quarter filling, the elec- 
tron density is one per two unit cells and the ground state 
is given by alternating electron-rich and electron-poor 
sublattice points with double degeneracy. If the spin de- 
grees of freedom of electrons are ignored for the moment, 
the electron-rich and electron-poor sites are mapped to 
spin up and down Ising variables; they are mapped to the 
Ising model. Therefore, at least at the classical level with- 
out quantum fluctuations, this free energy (1) provides a 
good starting point for general charge-order phenomena. 

In quantum systems with itinerancy of electrons, how- 
ever, a new aspect emerges, which is connected to the 
metal-insulator transition. In this paper, we study the 
first-order regime, where the metal-insulator transition 
occurs simultaneously with the charge-order transition 
as will be illustrated later in Fig. 17. We also study 
the continuous transition, where the charge-order transi- 
tion does not accompany the metal-insulator transition 
but only shows a crossover of transport properties. Then 
TCP appears also as the critical termination point of the 
metal-insulator transition. We will study the singularity 
of the conductivity in all the regions in detail. 

In this paper, we only consider common features of the 
charge-order transition with a commensurate periodicity 
and do not discuss details that depend on the detailed 
band structure, periodicity or charge-ordering pattern. 
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Criticality and tricriticality are well studied subjects of 
classical phase transitions, 5 while quantum effects par- 
ticularly on the tricriticality are not well explored. To 
study this issue using the microscopic model, we employ 
a simple extended Hubbard model defined by 



ij i 
(ij) ' 



(2) 



which captures some essence and offers useful starting 
point of studies for charge-order transitions, in general. 
Here, c ia (ci a ) creates (annihilates) an electron with spin 
a at site i, respectively and rii a and Ni = (n^ + n,jj are 
number operators. The transfer tij represents the elec- 
tron hopping between i and j sites and the terms pro- 
portional to U(V) represent onsite (intersite) Coulomb 
interactions, respectively. For the intersite repulsion V , 
we restrict ourselves to the repulsion of nearest-neighbor 
pairs for V. 

In the first step, we clarify phase diagram and critical- 
ity in the classical limit, where we take tij — 0. This clas- 
sical model partially captures essential aspects of con- 
tinuous charge-order transitions as well as TCP. After 
confirmation of the mean-field behavior, we will show re- 
sults of Monte Carlo studies by using the exchange Monte 
Carlo algorithm 9 to circumvent the critical slowing down. 
The critical exponents of various physical properties are 
studied. 

To consider quantum effects, we next study the itiner- 
ant model by taking nonzero tij by Hartree-Fock approx- 
imation and elucidate the phase diagram in the plane of 
the temperature T, the interaction couplings U and V 
and the parameter to control the band structure tij . We 
then clarify how the charge-order transition is coupled to 
the metal-insulator transition and how the criticality is 
described. 

The physical property and susceptibility studied in 
this paper is summarized in the following: 

1. Order parameter is defined as 



m = wX Nv 



(3) 



with N s being number of the sites and Q = (it, it) 
for the square lattice, for example. Order parameter 
susceptibility is defined as 



The susceptibility of to 2 is defined as 



(4) 



d 2 F 
d{m 2 ) 2 



2. Doublon density is defined as 

J_dF 1 ^ Ni(Nj-l) 
N s dU N s o 



(D(0)) 



and the doublon susceptibility is defined as 

d 2 F 
dU 2 ' 



Xd(o) 



(6) 



(7) 



The nearest-neighbor charge correlation which is con- 
jugate to the nearest-neighbor repulsion V is defined 

as 

<tj> 

and its susceptibility is defined as 

d 2 F 

XR 



dV 2 ' 

4. Internal energy is defined as 

F _ 1 d(F/T) 
N s d(l/T) 

and specific heat is defined as 

dE_ 
dT' 



(8) 



(9) 



C 



5. Charge density is defined as 

N s dn N s ^ 1 



n = — 



and charge susceptibility is defined as 



Xc = - 



d 2 F 



dfi 



2 ' 



(10) 



(11) 



(12) 



(13) 



The Mean-field values of critical exponents on the 
lambda line and TCP will be summarized in the next 
section in Table I. 

At TCP and on the lambda line, within the mean-field 
approximation, we will show that (D(0)), R and n have 
the same singularity as that of to 2 . Furthermore, Xd(o)> 
Xr and \c diverge with the same singularity as that of 
Xm 2 m general (see Table VI) . We will show that this re- 
lation still holds in Monte Carlo calculations. As shown 
below, only at half filling in the classical model, n does 
not have the same singularity as that of m 2 and there- 
fore Xc does not diverge because of the particle-hole sym- 
metry in both the mean-field approximation and Monte 
Carlo calculations. Strong divergence of various fluctu- 
ations at zero wavenumber is an outstanding feature of 
TCP in contrast to the lambda line. 

The organization of this paper is as follows: In §2, re- 
sults of mean-field and Monte Carlo studies of the clas- 
sical model are shown. We specify which types of fluc- 
tuations diverge at TCP in the classical model. In §3, 
we clarify effects of itinerancy on TCP and elucidate 
the criticality of the extended Hubbard model at the 
charge-order transition. In §4, we analyze the singularity 
of conductivity by Hartree-Fock approximation in detail 
by studying the singularity of the order parameter and 
compare the results with available experimental results. 
Section 5 is devoted to summary and discussion. 

2. Classical model 

2.1 Model 

First, we consider the classical limit of the extended 
Hubbard model. In electronic systems, Coulomb repul- 
sion plays a central role in stabilizing the charge or- 
der. Therefore this classical model, which only considers 
Coulomb repulsion, captures important aspects of the 
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charge order. We take tij = in eq. (2). Using the re- 
lation n^mi = (Ni(Ni — 1)) /2, we obtain the classical 
limit of the extended Hubbard model: 



(14) 



1 (ij) 1 

where Ni takes 0, 1 or 2. 

Here, we note that half filling n = 1 is realized by tak- 
ing /i = zV + U ~/2 from the symmetry of the Hamiltonian 
(14). Here, z is the coordination number. To see this, we 
make the transformation Ni — > 2 — Ni. This mapping 
changes the chemical potential /j, into p, = —fi + U + 2Vz, 
while at half filling this mapping does not change the 
total charge number. Therefore, we obtain the relation 
p, = fi. This leads to the relation fi — zV + U/2 for half 
filling. 

This model is equivalent to Blume-Emery-Griffith 
(BEG) model. 6 BEG model is used for analyzing the tri- 
critical behavior of 3 He- 4 He mixture, where the 3 He- 4 He 
mixture is realized in three-dimensional space. Since the 
upper critical dimension of TCP is three, the critical ex- 
ponents of TCP are correctly described by the mean-field 
treatment. The mean- field treatment indeed succeeded in 
explaining the phase separation of the 3 He- 4 He mixture 
and the existence of TCP. 6 

BEG model is defined as 



H. 



BEG 



(ij) » 



(15) 



where Si takes -1, or 1. The free energy of the BEG 
model is defined by 

^beg = ~T log Zbeg = -TlogTr[e- ffBEG / T ]. (16) 

Originally, the 3 Hc- 4 Hc mixture is realized in continuum 
space, while BEG model is its simplification to a lattice 
model, where 3 He and 4 He atoms occupy only discrete 
lattice points. Here, the Si = site corresponds to that 
occupied by 3 He atoms and Si = ±1 site corresponds to 
that occupied by 4 He atoms. Gauge degrees of freedom 
with 0(2) symmetry, which realizes superfludity of 4 He 
by its symmetry breaking, is drastically simplified to two 
discrete degrees of freedom Si = ±1. 

Using the correspondence relations Ni = Si + 1, U = 
A, V — — J and U/2+2V — \i = fi' , the classical extended 
Hubbard model (14) and BEG model (15) arc equivalent 
except for a constant. 

In BEG model, the order parameter is (Si). Here 
(Si) 7^ is interpreted as the superfludity of 4 He. In large 
A region, where 3 He concentration becomes higher, the 
superfluid phase becomes less favored. Therefore, in the 
large A region, the normal state can coexist with the su- 
perfluid state. Schematic phase diagram of normal and 
superfluid in the 3 He- 4 He mixture is illustrated in Fig. 3. 
Along the thick phase boundary, normal and superfluid 
phases coexist which means first-order transition across 
this phase boundary. The thin curve describes continu- 
ous transition line, which meets the first-order transition 
at the TCP. 

In BEG model, two different fluctuations diverge at 
TCP. One is the fluctuation of the order parameter. In 



T/J 



Lambda line 
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Superfluid f Normal Fluid 
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First-order transition line 



Fig. 3. Schematic phase diagram of BEG model. Thin line repre- 
sents the continuous transition line and thick line represents the 
first-order transition line. 



fact, x = — <9 2 Fbeg/<9// 2 diverges along the continuous 
transition and it also does at the TCP. The other fluctua- 
tion does not diverge along the continuous transition line 
but does only at TCP. This fluctuation is represented by 
the susceptibility of (Sf), as defined by 



X2 



d 2 n 



BEG 



<9A 2 



(17) 



The divergence of two fluctuations is an outstanding 
feature of TCP. In our classical model, A corresponds to 
U. Therefore the latter fluctuation corresponds to that 
of the doublon density (D(0)) = (n^riii), because n^riii 
is the conjugate quantity to U in (2). The divergence of 
the doublon susceptibility xd(o) defined in eq. (7) is in 
fact analyzed in the next section within the mean-field 
approximation. 

2.2 Mean-field approximation 

In this section, we show results of the mean-field ap- 
proximation just for comparison with Monte Carlo re- 
sults and results on the itinerant model in later sections. 
The mean-field approximation reproduces the Ginzburg- 
Landau approximation of the <f) 6 theory. In the mean-field 
approximation, we decouple the interaction term in eq. 
(14) as 

NiNj = (Ni - (Ni) + (Ni))(Nj - (Nj) + (Nj)) 

^N l (N 1 )+N 1 (N l )~(N l )(N J ). (18) 

We assume (A^) = n ± m. The upper sign + is for the 
A sublattice and the lower sign — is for the B sublatticc 
on a bipartite lattice. Using this mean field, we obtain 
the mean-field Hamiltonian 

H mf = j ]T Ni(Ni -l) + VzY / N l (n + m) 



i€A 



V z^Niin - to) - [i^Nt - 

ieB i 



N s Vz(n 2 - m 2 



(19) 

where z is the coordination number. From eq. (19), we 



J. Phys. Soc. Jpn. 



Full Paper 



Online-Journal Subcommittee 5 



obtain the free energy per site as 
I'M r = —TN s Zmf 



■— log(l + a(m) + 6(m)) 



-— log(l + a(-m) + 



V 



(n 2 -m 2 ), 



(20) 



where V = Vz and a(m) and 6(m) are defined as 

o(m) = e -(^(M^)-M)/r j (21) 



6(m) 



-(C7+2V'(n+m)-2/i)/T 



(22) 



The minimum condition dFMF/dm = leads to the 
self-consistent equation 

1 a(— m) + 26(— to) a(m) + 26(to) 
™ ~~ 2~4 + o(-m) + 6(-m) 1 + a(m) + 6(m) 

From this free energy, the charge density n and the 
doublon density (D(0)} are respectively given as 

9F MF 

1 a{m) + 2b{m) a(-m) + 2b(-m) 
2A +a(m) + b(m) 1 + a(-m) + b(-m) h K ' 

9Fmf 



n 



(D(0)) 



dU 



o(m) 



+ 



o(— to) 



■)■ 



2 v 1 + ci(to) + &(to) l + a(-m) + b(-m)' 

(25) 

We also obtain the internal energy E and the next 
nearest charge correlation R which is conjugate to V as 

d(F MF /T) 



E 



R 



d(l/T) 



U 



V'(n 2 +m 2 ) + ^(n-(D(0)))-tm, 

OFmf 
dV 

<»j> 

-(n -m ). 



(26) 



(27) 



For the moment, we consider half-filled case, n = 1. 
As mentioned above, half filling is realized by taking 
fj, = V' + U/2 independent of T or to. Therefore, we 
can easily expand the free energy with respect to the or- 
der parameter to at half filling. By taking n = 1 and 
/i = V + U/2, GL expansion is reproduced as 



F = A + rm 2 + um 4 + vm 6 - mH\ 



(28) 



where A, r, u and v are respectively defined as 

A = -V j2 - log(2 + e u / 2T ) (29) 
_ V'(2 + e u l 2T - IV IT) 



u = 



V = 



2{2 + e u / 2T ) 

V'\A-e u ' 2T ) 
12T 3 (2 + e c/ / 2T ) 2 

V 6 (-64 + 26e u / 2T -e u / T ) 



(30) 



(31) 



(32) 



360T 5 (2 + e u / 2T f 

and fft is the field conjugate to to. 
Lambda line: half filling 

Now, we consider the criticality of physical properties 
along the lambda line at half filling. As mentioned above, 
the lambda line is represented by the condition of r = 
and u > 0. Along the lambda line, the criticality of 
the order parameter and the susceptibility of the order 
parameter \m are defined by 

TO~r /3 , (33) 

Xm~r-\ (34) 

Since to = dF M F/dH^ and \ m — (d 2 F M f I '<9m 2 )~ t _ , 
we obtain from cq. (28) 



7=1. 



(35) 
(36) 



Along the lambda line, the criticality of (D(0)) in eq. 
(6) and Xd(o) m e Q- (7) are defined by 

(£>(0)} - r & 



(37) 
(38) 



In the mean-field scheme, the doublon density (D(0)) 
is proportional to to 2 except for a constant in the critical 
region. From eq. (25), we obtain the expansion of (D(0)) 
with respect to to as 



(D(0)) - a +fl!m 2 

where ao and oi are defined by 

1 

a = 



oi 



2 + e-^/ 2 ^ 

y/2 e J7/2T 



(39) 



(40) 



(41) 



2T 2 (l + 2e- l/ / 2T ) 2 ' 

Therefore, the singularity of the doublon density is 
equivalent to that of to 2 and its susceptibility Xd(o) is 
equivalent to Xm 2 - From cqs. (35) and (5), along the 
lambda line, critical exponents of the doublon density 
(-D(O)) and its susceptibility Xd(o) are obtained as 

02 = 1 (42) 

72 = 0, (43) 

where eq. (43) indicates that xd(o) does not diverge 
along the lambda line. 

Equations (26) and (27) show that E and R have the 
same singularity as that of to 2 whereas the specific heat 
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C and the susceptibility of R(xr) have the same sin- 
gularity as that of Xm 2 for the density fixed at n = 1. 
Therefore, along the lambda line, C and xr does not 
diverge. 

Here, we consider the criticality of the charge suscep- 
tibility Xc at the lambda transition. Equation (28) ex- 
presses the expansion when the density is fixed at n = 1 , 
while Xc is obtained by changing the density infinitcsi- 
mally from n = 1. Therefore, the singularity of Xc cannot 
straightforwardly be obtained from eq. (28). To clarify 
the singularity of Xc we calculate Xc numerically near 
the lambda line as shown in Fig. 4. We confirm that 
Xc does not diverge at the lambda transition within the 
mean-field approximation. 



0.44 



0.36 




0.28 0.3 0.32 

v/u 



Fig. 4. Charge susceptibility \ c is plotted near the lambda line 
at half filling. We fix temperature at T/U = 0.5. 

TCP: half filling 

Next, we consider the criticality of physical properties 
near TCP at half filling. TCP is represented by the con- 
dition of r = 0, u = in eq. (28). From this condition, the 
location of TCP is determined from eqs. (30) and (31) as 
T c = l/(41og2), V c = 3/(16 log2). At TCP, the critical 
exponents of the order parameter and the susceptibility 
of the order parameter Xm at TCP are defined by 

m ~ r A (44) 

Xm~r-T*. (45) 

Taking u = 0, we obtain from eq. (28) as 

b = \ (46) 
It = 1. (47) 

The critical exponents of (D(0)) and xd(o) at TCP 
are defined by 

(D(0)) ~ r &t (48) 
XD W ~r-^. (49) 

Again we note that, within the mean- field approxima- 
tion, the singularity of (D(0)) and xd(o) are also equiva- 
lent to that of m 2 and Xm 2 near TCP. Using (3 t = \ and 



the definition of x-m 2 m eq. (5), we obtain the critical 
exponents as 

fa = 2 $t = \ (50) 

72t = \. (51) 

From eqs. (26) and (27), the internal energy E and R 
have the same singularity as that of m 2 and their suscep- 
tibility have the same singularity as that of Xm 2 ■ From 
eq. (46), at TCP, we obtain the singularity of (D(0)>, E 
and R as 

(D{0)),E,R~ m 2 ~ r 02t , (52) 

X{D(o)), C, XR ~ Xm 2 ~ r~ l2t , (53) 

Equations (52) and (53) lead to (3 2 t = 1/2 and j 2 t = 1/2. 
Because j2t is positive, Xd(o), C and xr diverge at TCP 
in contrast to the lambda transition. 

Charge susceptibility at TCP 

We now consider the charge susceptibility Xc at TCP 
within the mean-field scheme. The charge susceptibility 
is defined by Xc = —d 2 FMF /d/i 2 . We confirm that Xc 
at half filling does not diverge at TCP because of the 
particle-hole symmetry (see section 2.4.3). Charge sus- 
ceptibility Xc at half filling is shown in Fig. 5. 

0.4624i 1 



0.462 
0.4616 
0.4612 




-20 -10 10 20 

(^xlO^ 

Fig. 5. Charge susceptibility % c is plotted near TCP at half filling 
within the mean-field approximation. We fix temperature at T c = 
l/(41og2). Critical interaction is given by V c = 3/(161og2). 

Next, we consider the TCP away from half filling. If 
we take fi/U = 1.2, we confirm that Xc diverges at TCP. 
Figures 6 and 7 show Xc at TCP. We estimate the crit- 
icality of Xc as Xc - \V - V c r ' 499(1) . This value is 
consistent with the mean- field value 72* = 1/2. This re- 
sult indicates that the singularity of the charge density 
n coincides with that of m 2 . This implies that the free 
energy can be expanded as 

F - A' - fi'n + v'n 2 + v'n 3 H (54) 

as a functional of n with A', //, u' and v' being constants. 
Therefore, Xc has the same singularity as that of Xm 2 at 
TCP, in general. 

Here, in Table I, we summarize the exponents of physi- 
cal properties obtained by the mean-field approximation. 
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Fig. 6. Charge susceptibility \ c is plotted near TCP at fi/U = 1.2 
within the mean-field approximation. Here, T c /U and V c /U are 
estimated as 0.31356 and 0.2663996, respectively. 
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Fig. 7. Scaling plot of Xc near TCP within the mean-field ap- 
proximation. Slope determines the exponent for \c- We estimate 
Xc = \V — V c \ ~ 0A 99W . This value is consistent with the mean- 
field value 724 = 1/2. 
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(3 2 = 2/3 = 1 


fo t = 2/3* = i 


E 
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Xm 
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7t = 1 
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at = 72 1 = 4 


XD(0) ~ Xm 2 


72 = 


1 i z — 

72t = 2 


XR 







Xc 


cusp 


cusp* 



Table I. Mean-field values of physical properties at half filling on 
the lambda line and at TCP. *Away from half filling (p,/U = 1.2), 
Xc diverges with the singularity X c ~ \ V - y c |-°- 499±0 - 001 . 



2.3 Monte Carlo results: Phase diagram 

To study this classical model beyond the mean-field 
approximation, we perform exchange Monte Carlo 9 cal- 
culations on the two-dimensional square lattice. We note 
that the Monte Carlo results give exact results within 
statistical errors. 

The exchange Monte Carlo calculation is performed in 
the following way: First, we prepare replicas that have 
different temperatures. On each replica, we perform a 
standard Monte Carlo simulation and renew each state. 
Then, we exchange the configurations of each replica by 
the extended probability distribution function W, which 
is defined as 



A(A„, X m ; T m ,T n ) = (l/T n -l/T n 



W(X n 



T 



1 



(E(X n ) - E{X m )) , 
(55) 

(56) 



1 + exp(A) ' 

where X n , E(X n ), T n represent the configuration of the 
nth-replica, the internal energy of the nth-replica, and 
the temperature of the nth-replica, respectively This ex- 
tended probability distribution function W satisfies the 
detailed balance condition 

,T n ) 



W(X n ,X m 



cxp(-A). 



(57) 



This method is efficient in overcoming the critical slowing 
down. 

At TCP, the upper critical dimension is at d c = 3. 
Therefore, the mean-field approach is justified in three 
dimensions except for logarithmic corrections. On the 
other hand, in 2D, the mean field theory is not accu- 
rate anymore. Here, we perform Monte Carlo calculation 
in 2 dimensions on a square lattice. 

Ground-state phase diagram is shown in Fig. 8. It 
should be noted that various charge-order states appear 
even in the classical model that considers only the onsite 
repulsion U and the nearest-neighbor repulsion V. 

At finite temperatures, we find TCP as in Fig. 9 at the 
phase boundary between the charge order with checker- 
board-type doublon alignment and the uniform phase. 
The order parameter of this charge order with checker- 
board-type doublon alignment at the wavenumber Q = 
(it, 7r) is given by 



' N 4 e 



(58) 



2.4 Monte Carlo results: Singularities and exponents 

We consider the case fj,/U = 1.2 as an example of a 
generic case. As shown below, \c at half filling has a 
different singularity from that at general fillings, while 
the singularities of other fluctuations are independent of 
fillings. Therefore, the case at fi/U = 1.2 captures the 
general nature of TCP. We will give the reason for the 
nontrivial behavior of Xc in detail later. 

First, we obtain the critical exponents v and 77 for 
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Fig. 8. Ground state phase diagram of the classical model. "Full 
filling" means a phase where all the sites are doubly occupied 
by particles. "Uniform" means a phase where each site is occu- 
pied by one particle. "3/4 filling" represents alternating singly 
and doubly occupied sites in a staggered pattern. Phase bound- 
aries between full filling and 3/4 filling, 3/4 filling and uni- 
form, uniform and quarter filling are given as fi/U = 8V/U + f , 
fi/U = 4V/U + f , fi/U = 4V/U, respectively. 



Exponent 
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2D Ising value 


V 
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2-du 
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? 


7 


(2 - r,)v 


4 



Table II. Hyperscaling relations on the lambda line, d represents 
spatial dimension. Exponents of 2D Ising model are obtained 
from exact solution. 11 



Exponents 


Hypcrscaling relations 


2D BEG model values 


vt 




0.56(f) 


Vt 




0.14(2) 


at 


2 - du t 


0.89(2) 


Pt 


2 


0.039(6) 


Pit 


v t d - 1 


0.12(2) 


It 


(2-r, t )vt 


1.03(2) 


fit 


2 — du t = a t 


0.89(2) 



Table III. Hyperscaling relations at TCP."' d represents spatial 
dimension. Exponents of 2D BEG model arc obtained numeri- 
cally. 10 The parentheses denote the error bars in the last digit. 



Lambda line 



0.11789 



Uniform 




Checker-board-tvpe 
doublort alignment 




0.2545 

First-order transition line 



On the lambda line, the universality class of this phase 
transition is categorized to that of the conventional 2D 
Ising model. 7 As described below, we indeed confirmed 
that critical exponents are the same as that of the 2D 
Ising model derived from Onsager's exact solution. 11 

To obtain critical exponents, we use finite-size 
scaling. From eqs. (59) and (60), at t = 0, 
log[M(Li)/M(L 2 )]/ log[ii/L 2 ] has the universal value 
2 — r/, which is independent of lattice size. Similarly, 
\og[A(Li) / A(L2)]/ log[Li/ L 2 ] converges to its universal 
value \ jv. 

The result of finite-size scaling at V/U = 0.6 is shown 
in Figs. 10 and 11. We obtain rj = 0.25(1) and v — 
1.03(5). We estimate the error by the standard devia- 
tion of the crossing points. These results are consistent 
with the exact values r\ = 1/4 and v = 1. 



Fig. 9. Phase diagram in the parameter space of scaled temper- 
ature T/U and the intersite interaction V/U at fi/U = 1.2. 



finite-size scalings defined by 

M — (to 2 ) = L 2 ^ v ^i(tL 1 ^), 



d\n{m 2 } 



(59) 
(60) 



Here, L is the linear dimension of the lattice size, t is 
the reduced temperature defined by t = (T — T c )/T c and 
$i and $2 are scaling functions. As shown in Table II 
and III, using hypcrscaling relations, 5 other exponents 
defined in Table I are described by using v and to 

2.4-1 Exponents v and r\ analyzed from finite-size scal- 
ing 

Lambda line: fi/U = 1.2 



I 

130 
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£,= I28,l r l6 ■ 
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■ I[= 128, 1 2 =32- 
l t = 128, 1 2 =48- 

■ii= I28.£ 2 =<54- 



0.6 



0.604 



0,608 



T/U 



Fig. 10. Result of finite-size scaling near the lambda line. The 
parameters arc L\ = 128, L2 = 16 (open diamonds), L\ = 
128, L 2 = 24(+), Li = 128, L 2 = 32 (squares), Li = 128, L 2 = 
48(x) and L\ = 128,1/2 = 64 (triangles). The horizontal line 
shows the exact critical exponent. 
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0.1179 
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0.118 



Fig. 11. Result of finite-size scaling near the lambda line. The 
parameters are L\ = 128,1/2 = 16 (open diamonds), L\ = 
128, L 2 = 24(+), Li = 128, L 2 = 32 (squares), L 1 = 128, L 2 = 
48(x) and L\ = 128,1,2 = 64 (triangles). The horizontal line 
shows the exact critical exponent. 



Fig. 13. Finite-size scaling of A near TCP. The parameters are 
Li = 128, L 2 = 16 (open diamonds), L 1 = 128, L 2 = 32(+), 
Ll = 128, L 2 = 48 (squares), Li = 128, L 2 = 64(x). The hori- 
zontal line shows the critical exponent in Rcf. 10. 



TCP: li/U = 1.2 

We next consider the criticalities of TCP. We perform 
finite-size scaling for M and A. We determine the loca- 
tion of TCP by the best fit of these finite-size scaling. 
We obtain the best fit for M and A at V/U = 0.2545. 
By the finite-size scaling of M (Fig. 12), we estimate 
the critical temperature T c /U = 0.111790(1) and the 
critical exponent rj t = 0.11(2). By the finite-size scal- 
ing of A (Fig. 13), we estimate T c /U = 0.111788(1) 
and the critical exponent v t = 0.55(3). These tricriti- 
cal exponents are consistent with those obtained by the 
Monte Carlo renormalization group (MCRG) method. 10 
Results by the MCRG method in the literature 10 shows 
rit = 0.14(2), v t =0.56(1). 




1178 



1179 

T/U 
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Fig. 12. Finite-size scaling of M near TCP. The parameters are 
Li = 128, L 2 = 16 (open diamonds), Li = 128, L 2 = 32(+), 
Ll = 128, L 2 = 48 (squares), L\ = 128, L 2 = 64(x). The hori- 
zontal line shows the critical exponent in Ref. 10. 



We summarize the critical exponents v and r\ obtained 
by Monte Carlo calculations in Table IV. 

2.4-2 Various singularities and exponents 

We now clarify the nature of fluctuations that diverge 
at TCP and the lambda line. In this section, we will 
consider five fluctuations: Xm, Xd(o)j XR an d Xc- 
At TCP, the critical exponents j t and j2t are defined 





V 


V 


Lambda line 


0.25(1) 


1.03(5) 


TCP 


0.55(3) 


0.11(2) 



Table IV. Critical exponents v and rj on the lambda line and at 
TCP. The values in parentheses denote the error bars in the last 
digit. 



by 



X. 71 



Xd(o) 



\T- 
\T- 



T c 



(61) 
(62) 



Using the hyperscaling relation in Table III, we estimate 
7t and j 2 t as 



7 t = (2 - nt)v t = 1-04(6) 
72* = 2 - 2v t = 0.90(6). 



(63) 
(64) 



In Monte Carlo calculations, we will show that other fluc- 
tuations C, xr and Xc has the same singularity as that of 
Xd(o) m general. We consider Xc at half filling separately 
in the next section. 

Lambda line: li/U = 1.2 

On the lambda line, the universality class of the phase 
transition is described by that of the 2D Ising model. 
Therefore, Xm behaves as Xm ~ \T — T c |~ 7 ,7 = 7/4, and 
the singularity of the specific heat is logarithmic, C ~ 
log \T — T c \. This logarithmic singularity of the specific 
heat causes logarithmic singularities of other fluctuations 
such as Xd(o) and xr as proven by Ehrenfest law (see 
Appendix A). From cqs. (89), (90), (91) and (92), Xd{o), 
Xr and Xc diverge as log \T - T c \. 

TCP: fi/U = 1.2 

From the hyperscaling relation in Table III, it is clear 
that the specific heat C has the same critical exponent 
as that of the doublon density susceptibility Xd(o) - How- 
ever, the singularities of Xc and xr are still not clear. 

To apply Ehrenfest law, it is necessary to define partial 
derivatives. At the endpoint of the continuous transition 
line, for example at TCP, partial derivatives are not well 
defined. Therefore, at TCP, Ehrenfest law does not apply 
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T c 


exponent 


Xd(o) 


0.11786(3) 


1.7(2) 


Xr 


0.11786(3) 


1.7(2) 


C 


0.11785(3) 


1.5(2) 


Xc 


0.11784(2) 


1.4(2) 


72t 
ft 




1.6(1) 



Table V. Critical temperatures and critical exponents at TCP 
obtained from Monte Carlo calculations. The values in paren- 
theses denote the error bars in the last digit. 

in general. From this, it is nontrivial whether xr an d Xc 
have the same singularity as that of Xd(o)- 

Near TCP, we confirm that the doublon susceptibility 
Xr>(o)> XRi the specific heat C, and the charge suscepti- 
bility Xc are enhanced. Using finite-size scaling for Xn(o): 
Xr, C and \c, we estimate critical temperatures and the 
critical exponent ^t/vt as shown in Table V. These re- 
sults are consistent with the value estimated from the 
hyperscaling relation in Table III as 

l2t /v t = 2/v t - 2 = 1.6 ± 0.1. (65) 

From these results, we conclude that \R and Xc have 
the same singularity as that of Xd(o)- The positive 72* 
indicates that Xd, XR> C and Xc diverge at TCP. 

2.4-.S Charge susceptibility: half filling 
Lambda line 

At half filling, Xc is not enhanced near the lambda 
line, while other fluctuations are enhanced. This is due to 
the particle-hole symmetry. Because of the particle-hole 
symmetry, the phase diagram should be symmetric in the 
/z-T plane, as shown in Fig. 14. Half filling corresponds to 
the peak of the lambda line and the slope of the lambda 
line at half filling is zero. From the condition that the 
slope of the lambda line in the \x-T plane is zero and eqs. 
(87) and (88), Xc has a singularity weaker than that of 
the specific heat C leading to a singularity weaker than 
the logarithmic one on the lambda line, and we do not 
observe the enhancement of Xc on the lambda line. 



T 








Fig. 14. Phase diagram in the fi-T plane for V/U > 0.2546. 
Hhal f = 4V + U/2 is the chemical potential that corresponds 
to the case of half filling. The transition point at half filling is at 
the peak of the lambda line in the /i-T plane. 



TCP 

At half filling, we estimate the location of TCP at 
T c /U = 0.1567 and V c /U = 0.2546 by the best fit for 
M and A in eqs. (59) and (60). 

As mentioned above, Ehrenfest law docs not apply to 
TCP in general. However, only at half filling Ehrenfest 
law can be applied to TCP, because TCP is not the end- 
point of the lambda line in the [i-T plane as shown in 
Fig. 15. Using the same reasoning for the lambda line, 
Xc has a singularity weaker than that of the specific heat 
C. Therefore Xc does not diverge with the singularity 
\T — T c \~ l2t and we do not find the divergence of Xc near 
TCP at half filling, as shown in Fig. 16. 



T 




► 



Fig. 15. Phase diagram in the /i-T plane at V/U = 0.2546, which 
contains TCP. TCP at half filling is located at the peak of the 
lambda line in the fi-T plane. 



0.5 




0.2 1 1 

0.15 0.16 0.17 

T/U 

Fig. 16. Xc near TCP (V c /U = 0.2546, T c /U = 0.1567) at half 
filling. The parameters are V/U = 0.2546 (open diamonds), 
V/U = 0.2548(+), V/U = 0.2550 (squares) and the linear dimen- 
sion of the lattice size is L = 128. We do not find the divergence 
of Xc up to L = 128. 

We have clarified which fluctuations diverge at TCP in 
both mean-field approximation and Monte Carlo calcula- 
tion. We confirm that such divergences occur in the same 
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MF TCP 
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/3 = 0.12(6) 
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h'Zt 9 




Sot = 10(61 


Xm 


7 = 1 


7t = 1 


7 = 1.8(1) 


7t = 1.04(4) 


Xd(o) 


72 = 


72t = k 


log 


72t = 0.9(2) 
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a = 




log 


at = 0.8(2) 


XR 







log 


0.9(2) 


Xc 


cusp 


cusp* 


log** 


0.8(2)** 



Table VI. Critical exponents at lambda line and TCP. MF and 
MC represent mean-field results and Monte Carlo results, respec- 
tively. The mean-field results are obtained at half filling and the 
Monte Carlo results are obtained in the case at fi/U = 1.2. The 
values in parentheses denote the error bars in the last digit. * 
Away from half filling (fJ./U = 1.2), \c diverges with the singu- 
larity xc ~ I V - V c |-°- 499 . "At half filling, X c docs not diverge 
on the lambda line and at TCP. 



quantities between the mean-field approximation and the 
Monte Carlo results. Exact calculations give only modifi- 
cations of the critical exponents from those of the mean- 
field results. Except for the charge susceptibility Xc ah 
the fluctuations we consider diverge with the power law 
at TCP. Because of the particle-hole symmetry, \ c does 
not happen to diverge at half filling. However, away from 
half filling, we confirm that Xc generically diverges at 
TCP with the same singularity as that of Xd(o)- Finally, 
the critical exponents we obtain in the classical model is 
shown in Table VI. We note that various susceptibilities 
and the specific heat diverge much more strongly than 
those at the lambda transition. We have for the first time 
shown that the critical exponent of the charge suscepti- 
bility Xc and the susceptibility of the next-nearest cor- 
relation xr is the same as that of the susceptibility of 
Xd(o) m general. It is remarkable that some fluctuations 
at zero wavenumber diverge as Xd(o) m contrast to the 
lambda line. 

3. Itinerant models 

3.1 Model 

In this section, we consider the effect of itinerancy of 
particles on TCP within the Hartree-Fock approxima- 
tion. We add the hopping term to the classical model. 
In this section, we study charge ordering transition as 
well as the competition between metals and insulators. 
Metal-insulator transitions are also discussed in terms of 
the singularity of conductivity. We first consider the ex- 
tended Hubbard model with the next-nearest-neighbor 
hopping(t') on the two-dimensional square lattice. The 
Hamiltonian is given by eq. (2). 

3.2 Hartree-Fock approximation 

We consider checker-board-type doublon alignment de- 
scribed by the order parameter to defined from m = 
1 / N sJ2i(( n il + n J i)cxp(iQr i )). Here, we take (n iT ) = 
i n il) = (n + TOCxp(zQri))/2, with n being the average 
charge density given by l/iV s X)i( n iT + n ii)- Within the 
Hartree-Fock approximation, we decouple the interaction 
term. 



The onsite interaction term is decoupled as 

i 

i 



U 



E 



Uia ( n + me 



iQi-i 



)) 



UN, (n 2 + to 2 ) 



(67) 



The nearest-neighbor interaction term is decoupled as 

^E^i ( 68 ) 

m 

= Vj2m (N t ) + (Ni))(Nj (N 3 ) + (Nj))) 

«> 

= AnV^N, - AmV ^ N t exp(iQri) - 2VN s (n 2 - to 2 ). 

i i 

(69) 

Finally, we obtain a Hartree-Fock Hamiltonian (Hhf) 
in momentum space as 

Hhf = E + it + 4nV ) c ^ Ckff (70) 

_ ^N s{ n^rn 2 ) +WNs{n2 _ m2) y (n) 

where e(k) = —2t(cos(k x ) + cos(k y )) + At' cos(k x ) cos(k y ) . 

Diagonalizing the Hamiltonian leads to two bands of 
the form 

e(k)+e(k + Q) Un 



E± = 



+ — +AnV (72) 



*v( c(k) -f +Q) )* + W, 

where g = AV — -j . 

Using this Hartree-Fock band dispersion, we obtain the 
free energy 

T 

Fhf = -—\ogZ HF 



.^5>g(l + e-<*. 

k,r/ 

U(n 2 + to 2 ) 



(k)- M )/TN 



+ 2V(n 2 - to ) 



(73) 



where the suffix r\ takes ±. 
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From this free energy, the charge density n and 
the order parameter m are determined from the self- 
consistency condition 

2 



— V/(£ + (k)) + /(E_(k)) 



k 



(74) 



/(£L(k))-/(£ + (k)) 



£+(k) 



£7_(k) 



(75) 



where /(a:) is the Fermi-Dirac distribution function. 

We also obtain the doublon density (D(0)) and R con- 
jugate to V defined by 

(£>(0)) = <n iT >(n, 1 > = ^±^ ) 



R = (NJiNj) = 2{n 2 - to 2 ). (76) 

(D(0)) and R are described by the square of to, which 
are notable features of the mean-field approximation. 

It turns out that the self-consistent equation (75) with 
cq. (73) for the extended Hubbard model (2) is equiva- 
lent to the self-consistent equation for the simple Hub- 
bard model obtained by taking V — in eq. (2). The 
equivalence holds by replacing the charge-order param- 
eter m with the antiferromagnetic(AF) order thaf — 
Ei(( n »T _ "-»i) ex P(Q r i))/^ and by putting g = U/2. 
However, one should be careful about this mapping. This 
mapping is complete only at the Hartree-Fock level. The 
charge order is the consequence of discrete symmetry 
breaking while the AF order is realized by the contin- 
uous symmetry breaking of SU(2) symmetry. Therefore, 
if it would be exactly solved in two-dimensional systems, 
the charge order can indeed exist at finite temperatures, 
though AF order cannot exist at finite temperatures by 
Mermin-Wargner theorem. Therefore, the Hartree-Fock 
approximation captures some essence of the charge order 
at finite temperatures, while the AF order at nonzero 
temperatures in two dimensions is an artifact of the 
Hartree-Fock approximation. 

3.3 Critical exponents 

We first consider the case with the electron density 
fixed at n = 1 in the canonical ensemble. We determine 
the location of TCP, and the critical exponent of the 
order parameter. The behavior of the order parameter 
m at TCP is different from that of the lambda line as 



m oc \g-g c \ 



13 



h for lambda line 



Pt = \ for TCP 



A schematic phase diagram is shown in Fig. 17. The 
TCP line appears at finite temperatures in T-g-t' space. 
In this paper, we study TCP at t' jt = 0.2 in detail. 
The phase diagram at t' jt — 0.2 is shown in Fig. 18. The 
singularity of the order parameter is shown in Fig. 19 and 
scaling result is shown in Fig. 20. The critical exponent of 
the order parameter fit = 0.249 ± 0.002 is well consistent 
with the expected mean-field value (3t = 1/4. 

We specify which fluctuations diverge at TCP and on 
the lambda line in the itinerant model (2). The singular- 
ity of the doublon susceptibility (Xd(o)) an d fluctuation 
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Charge ordered phase 




Fig. 17. Phase diagram in T-g-t' space. The lambda (first-order) 
surface is described as green (red) surface above (below) solid 
and broken lines, respectively. The TCP line is shown in blue. 
Near zero temperature, at the yellow circle, the TCP line goes 
into the ordered phase and changes into the critical line of the 
metal-insulator transition (broken black line). 



0.06 



0.04 



0.02 




First-order transition line 



Fig. 18. Phase diagram at t'/t = 0.2. We estimate the location 
of TCP as Tc/t = 0.0511 and g c /t = 0.01030 



conjugate to the nearest-neighbor repulsion V (xr) are 
respectively given as 



Xd(o) 



d(D(0)) 
dU 



1.9 - 9c 



2/3-1 



OR . 1213—1 

At TCP, since (3 t = 1/4, Xd(o) an d Xr respectively be- 
have as 



Xd(o) ~ \9-9c\ 2 

XR ~ 1.9 ~ ffc| - '- 

This means that Xd(o) an d XR diverge at TCP, in con- 
trast to the criticality on the lambda line. 

The singularity of the charge susceptibility will be dis- 
cussed in §3.6. 

3.4 Metal-insulator transition 

We consider the relation between metal- insulator tran- 
sitions and charge-order transitions at t'/t — 0.2 in 
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10 4 x(g-g,)/r 

Fig. 19. Singularity of the order parameter at TCP. A is de- 
fined by A = mg, which has the same singularity as that of 
the order parameter m. We estimate tricritical temperature as 
T c /t = 0.0511 and tricritical interaction as g c /t = 0.01030 



6 




io 4 *(g-&)/f 

Fig. 20. Critical exponent of the order parameter at TCP. A has 
the singularity such that A = a\g — <; c |' 3t • Therefore we obtain 
(dlog A/dg) -1 = \g — g c |. Using this relation we estimate 
fit = 0.249 ± 0.002 as the fitting shown by the solid line. 

this section. At zero temperature, the metal-insulator 
transition and the charge-order transition occur at the 
same time. Namely, the charge gap opens for the 
whole Brillouin-zone simultaneously with the charge- 
order transition. The critical A of the metal-insulator 
transition at t'/t = 0.2 is given by A c /t = 2t'/t = 0.4. 
The gap A shows a jump from A/t = to A/t = 0.4031 
at the first-order transition point at g/t — 1.0319. Figure 
21 shows the band dispersions of the parametal side and 
charge-ordered side at the transition point. 

At finite temperatures, the jump of A decreases and 
vanishes at TCP (T c /t = 0.0511). Therefore, the metal- 
insulator transition, strictly speaking metal-semimetal 
transition, and the charge-order transition occur sepa- 
rately at high temperatures. Figure 22 shows the band 
dispersions in the parametal side, -Ef ree and in the charge- 
ordered side, E± at the transition point at T/t = 0.04. 
The gap A shows a jump from A/t = to A/t = 0.2170 
at the first-order transition point at g/t = 1.0310. There- 
fore, the charge gap does not fully open for the whole 
Brillouin zone. 
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Fig. 21. Band dispersions at the transition point at zero tem- 
perature are shown. The Fermi energy for the charge-ordered 
side is described by fip. Band dispersions are given by 
E± = At' cos(k x ) cos(ky) ± y/4t 2 (cos(k x ) + cos(fc y )) 2 + A~| for 
the charge-ordered phase and Ef ICC = —2t(cos(k x ) + cos(fc y )) + 
At' cos(fc^) cos(ky) for the parametal phase where A c /t = 0.4031. 
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Fig. 22. Band dispersions at the transition point at T/t = 
0.04 are shown. The Fermi energy for the charge-ordered side 
is described by fip. Band dispersions are given by E± = 
4t' cos(fc^) cos(fc a )± i/4t 2 (cos(fc I ) +cos(fc !/ )) 2 + A 2 and i?f rco = 
-2t(cos(fe a ,)+cos(fcj / ))+4t' cos(fc;c) cos(fc !/ ) where A c /t = 0.2170. 

3. 5 Conductivity 

In this section, we consider the singularity of carrier 
density. The definition of carrier density is given by the 
electron density in the upper band as 

X= [ f(E+(K))dk. (77) 
Jbz 

Phenomenologically the conductivity may be expressed 
by the product of the carrier density X and the carrier 
relaxation time r as a cx Xt. In this section, we consider 
the singularity of the carrier density X. The singularity 
of the relaxation time will be considered later. 

First, we consider how the carrier density depends on 
the gap A. In our calculation, we first fix A and deter- 
mine the interaction g using the self-consistent equation. 
This allows us to determine the relation between the car- 
rier density and A even in the first-order transition re- 
gion. From eq. (101) in Appendix B, the singularity of 
carrier density is given as \X — X c \ = A 2 (^4 log A + B) 
in the asymptotic region. This relation holds both near 
the lambda transition as well as near TCP. 
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On the lambda line, the conductivity exponent p de- 
fined by \X - X c \ oc (g - g c ) p (A A log A + B A ) is given 
from p = 2(3 = 1 since A oc (g — gc)' 3 with (3 = 1/2. 

At TCP in the canonical ensemble, the conductivity 
exponent p t defined by \X — X c \ oc (g — g c ) Pt (A\og A+B) 
is given from p t = 2(3 t = 1/2 since A oc (g — 5 C ) /3 * 
with (3 t = 1/4. As we will show later, at TCP in the 
grand-canonical ensemble, the critical exponent (3 t = 1/4 
changes into (3 = 1/2. Therefore the conductivity expo- 
nent pt changes into p' t = 2/3 = 1. The exponent p on 
the lambda line is twice as large as p t at TCP in the 
canonical ensemble. 

Numerically, near TCP, we confirm that carrier density 
has the singularity defined in eq. (101). The singularity 
of the carrier density is shown in Fig. 23. We fit (X c — 
X)/A 2 with the function Alog A + B, then obtain A = 
-0.276(2) and B = -0.73(2). This result is shown in Fig. 
24. 
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Fig. 23. Carrier density X is plotted near TCP. 
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Fig. 24. Singularity of (X c - X)/A 2 is plotted near TCP 
(T/t=0.0511) as (red) diamonds. We confirm that (X c -X)/A 2 
has the singularity A log A + B shown by the solid line, where 
we estimate A = -0.276(2) and B = -0.73(2). 
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3.6 Phase Separation 

In this section, we consider the phase separation near 
TCP in the grand-canonical ensemble. In our model, we 
confirm that the phase separation occurs at U/t = 10.0 
not only along the first-order transition line but also 
along the second-order transition line as shown in Fig. 
25. Therefore TCP in the canonical ensemble is actu- 
ally unstable toward the phase separation if the grand- 
canonical ensemble is employed. If we calculate in the 
grand-canonical ensemble, the location of TCP shifts 
from the canonical one. 

Here, we note that the region of the phase separation 
depends on the onsite interaction U. The definition of /j, 
is given as 

M = MF + (4V+^)n 

= HF + (g + U)n (78) 

where /if is the Fermi energy and it only depends on tem- 
perature and band structure, namely t, t' and A. Since 
U is independent of g, fi depends on U and the region 
of the phase separation depends on U. Roughly speak- 
ing, the phase separation region shrinks for larger U and 
above the possible threshold U c , TCP of the canonical 
ensemble becomes stable. 

At U/t = 10, the location of TCP is at T c /t = 0.05726 
and g c /t = 1.02905. Because the half-filling plane is tan- 
gential to the lambda surface in T-fi-g space near TCP 
as we see in Fig. 26, the critical exponent of the order 
parameter is different from the generic one of TCP but 
the same as that of the lambda line. In the route that 
approaches TCP along the lambda line, the critical ex- 
ponents of the order parameter and the charge suscepti- 
bility are respectively given by 

m~\ g - gc f (79) 

Xc~\g-9c\-^ (80) 

where (3=1/2 and 72 = l. 5 However, for the generic 
filling-control transition toward TCP, it is away from the 
lambda line and the exponents should recover the generic 
tricritical exponents (3 t = 1/4 and 72* = 1/2. 

Figure 26 is the phase diagram at the critical temper- 
ature T c /t = 0.05726. We obtain the singularity of the 
order parameter as m ~ \g — g c | - 499 ± - 002 This value is 
consistent with the mean-field value. However, the expo- 
nent pa is the same as that of the generic one. There- 
fore, we obtain the singularity of the carrier density as 
\X - X c \ ~ \g - g c \ p 't where p' t = 2 x /3 = 0.92 ± 0.01. 
In the present model, TCP appears to move out of the 
phase separation region for U/t > 100, and canonical and 
grand-canonical results give the same TCP. 

In real materials, long-range Coulomb interaction 
plays an important role and it can suppress the phase 
separation. Therefore, the TCP in real materials may be 
compromised and located between the grand-canonical 
and canonical TCP's. Critical exponents in this case 
are not clear at the moment. Effects of the long-range 
Coulomb interaction on the phase separation is left for 
future studies. 
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Fig. 25. Phase separation region in the case of U/t = 10 is de- 
scribed as pink hatched region. TCP in the canonical ensemble 
is unstable and TCP in the grand-canonical ensemble shifts to 
the point at T c /t = 0.05726, g c /t = 1.029046. 
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Fig. 26. The lambda line and the line which corresponds to the 
half filling in the fi-g plane at the critical temperature T c /t = 
0.05726 are plotted. Near TCP, half-filling line becomes tangen- 
tial to the lambda line. Therefore, the singularity of the order 
parameter is given by m ~ \g — g c \@ where j3 = 1/2. 



We calculate the charge susceptibility Xc near TCP 
in the grand-canonical ensemble. Figure 27 shows the 
charge susceptibility Xc near TCP at T c /t = 0.05726. 
We obtain the singularity of charge susceptibility as \c ~ 
\g — g c \~ 12 where 72 = 0.99 ± 0.01 as shown in Fig. 28. 
This value is consistent with the mean-field value 72 = 
l. 5 

4. Comparison with experimental results 

In this section, we compare the present Hartree-Fock 
analysis with experimental results. An organic conduc- 
tor (DI — DCNQI)2Ag is a compound with the uniform 
stacking of planar DI-DCNQI molecules along the c axis 
forming conductive chains with monovalent and non- 
magnetic Ag ions stacked among the chains. This com- 
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Fig. 27. Xc near TCP in the grand-canonical ensemble is plotted. 
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Fig. 28. Scaling analysis of \ c . Slope of green line corresponds 
to the critical exponent 72. We estimate 72 = 0.99 ± 0.02. This 
value is consistent with the mean- field value 72 = 1. 



pound has a quasi-one-dimensional quarter-filled 7r-band 
system. (DI — DCNQI^Ag undergoes a transition to a 
charge-ordered insulator, where the 3:1 charge dispropor- 
tionation occurs on alternate molecules within a chain. 
This compound shows continuous and first-order phase 
boundaries together with TCP on the phase boundary 
between metals and charge order drawn in the phase di- 
agram in the plane of temperature and pressure. The 
pressure is supposed to control the bandwidth t in our 
model and the extended Hubbard model (2) offers a min- 
imal model and may be a good starting point to discuss 
the qualitative aspects of the phase diagram as well as 
the critical behavior. The three dimensional charge order 
takes place in this compound. The upper critical dimen- 
sion of TCP is three, so it is expected that the singularity 
of conductivity is well described by that of the Hartree- 
Fock approximation. The universality and critical expo- 
nents at charge-order transitions should not depend on 
details of the charge-order structure and the estimated 
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exponents in the extended Hubbard model in the previ- 
ous sections must apply here. 

Resistivity has been measured near TCP by Itou et al. 4 
and they have estimated that the location of TCP is near 
18 kbar. We analyze the singularity of conductivity at 18 
kbar and 17 kbar. Experimental data of the conductivity 
is shown in Fig. 29. 
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m) 

Fig. 29. Experimental data of the conductivity near TCP. 4 



singularity \<t-<j c \ oc |T-T c |«(Alog \T-T C \ + B), where 
q = p t = 1/2 (g = p' t = 1) in the canonical (grand- 
canonical) ensemble, respectively. To estimate q, we plot 

"{^ a^tUi + b ) '■ < 81 > 

We note that the slope of S becomes q~ x in the region 
where \T — T c \ <C 1. The results of scaling analysis are 
shown in Fig. 31. 
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We differentiate the conductivity with respect to tem- 
perature T. If the conductivity has a singularity such 
as a - a c - \T - T c \ q (Alog\T -T c \ + B), the dif- 
ferentiation of the conductivity behaves as da/dT oc 
\T - T^- 1 log |T - T c \. If q is less than 1, da/dT should 
diverge at T c . Even if q is equal to 1, da/dT has loga- 
rithmic divergence. Data of da/dT is shown in Fig. 30. 
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Fig. 30. Experimental data of da/dT near TCP. 4 



We assign the critical temperature T c as the mid- 
dle point between the highest and the second highest 
points of da/dT. Using this T c , we determine a c , and 
perform scaling analysis. We estimate T c — 28.23 K, 
a c = 1260.13ohm" 1 cm" 1 for 17 kbar and T c = 23.49 
K, a c = 992.37ohm" 1 cm" 1 for 18 kbar. 

From the analysis of Hartree-Fock approximation, if 
the singularity of conductivity comes from that of the 
carrier density, it is expected that conductivity has the 



Fig. 31. Scaling analysis of a at 17 kbar and 18kbar in the metal- 
lic and insulating sides. The slope of the broken straight line 
without symbols corresponds to the critical exponent of conduc- 
tivity q = p' t = 1 without logarithmic correction. We estimate 
slope as 1.03(6) for 17 kbar and 1.1(1). 

In both data for 17 and 18 kbar, slopes become 1 
near T c . These results indicate that the phase separa- 
tion occurs along the first-order transition line and ter- 
minates at TCP in (DI - DCNQI) 2 Ag. In this case, the 
Hartree-Fock analysis in the grand-canonical ensemble is 
justified in three dimensions and the critical exponent 
is q = p't = 1. Therefore the experimental data for 17 
and 18 kbar are sufficiently close to the critical region 
of TCP in the grand-canonical ensemble. We do not ex- 
clude another possibility that the phase separation does 
not occur or terminates on the first-order transition line 
of the canonical ensemble in (DI — DCNQI^Ag. In this 
case, Hartree-Fock analysis in the canonical ensemble is 
justified in three dimensions and the critical exponent 
is q = p t = 1/2. Therefore, the experimental data are 
not sufficiently close to the critical region of TCP. If one 
approaches the critical region of TCP, the singularity of 
the conductivity should approach q = p t = 1/2 if the 
phase separation is absent. More complete understand- 
ing is left for future studies. It is desired to measure the 
critical properties of the conductivity in the region closer 
to TCP after identifying TCP precisely. It is also intrigu- 
ing to observe the degrees of carrier-density inhomogenc- 
ity in samples, which is in fact helpful in revealing the 
tendency for the phase separation. 

Here, we consider the singularity of relaxation time. 
Phenomelogically, carrier density X (relaxation time r) 
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is divided into the coherent part X co h{T co h) and the in- 
coherent part X inc (T inc ). Then, the conductivity may be 
expressed as a ~ X co hT co h + X inc T inc . If we consider 
the singular part of X coh (X inc ) and T coh (T inc ) associ- 
ated with the transition as X s coh {X s inc ) and t c s o/i (t/„ c ) , 
such that X coh = X° coh + X s coh (X mc = X? nc + X? nc ), 
Tcoh = t c ° o/i + T s coh {rl nc = r°„ c + r?„ c ) with the normal 
parts X° oh (X? nc ) and T c ° oh (T°„ c ) , we obtain X coh r coh ~ 

X° T S -4- X S T° 4- X S T S (X- T- ~ X° T S 4- 

x !nc T inc +^i f iic T Di respectively. First, we consider the 
singularity of the coherent part. It is expected that X^ oh 
has the same singularity as that of X. Namely, the sin- 
gularity of X s coh is given as \T-T c \*(Alog \T - T c \+B). 
In Appendix C, it is speculated that r co h either does not 
have the normal part T® oh or does have the nonzero T® oh . 
Here, r^ oh has the same singularity as that of the carrier 
density X except for logarithmic correction. From this, 
the dominant singularity of X^^t^ comes from r s coh or 
X s coh and the singularity of T* oh is weaker than that of the 
carrier density. Next, we consider the singularity of the 
incoherent part. Again, it is expected that X? nc has the 
same singularity as that of X. One of the evidences that 
rf nc has a weaker singularity than that of the carrier den- 
sity is found in the study of the double-exchange model 
which describes ferromagnetic transition in perovskite- 
type manganese oxides. Theoretically, the singularity of 
the relaxation time near the ferromagnetic transition is 
studied in the double-exchange model within the mean- 
field framework. 15, 16 In this system, carrier density does 
not change at the ferromagnetic transition, and the sin- 
gularity of the conductivity comes from that of the re- 
laxation time. The singularity of the relaxation time T? nc 
is characterized as that of the square of the order pa- 
rameter m (spontaneous magnetization). Therefore, r? nc 
is given as 

T- cK m 2 cx|T-T c r. (82) 

Experimentally, for example in Lai^Sr^MnOa, the sin- 
gularity of the resistivity is also characterized as that of 
m 2 near the transition point. 17 This singularity is weaker 
than that of the carrier density given in cq. (101). This 
argument is for the case of ferromagnetic transition and 
it is not clear how it is universal. However it is likely 
that in the case of the charge-order system it shares the 
similar mechanism and r? nc has the similar singularity. 

We finally obtain the singularity of a as the same one 
with X. Detail analysis of the singularity of the relax- 
ation time in the charge-order system is left for future 
studies. 

5. Summary and Discussion 

A simple charge-order transition has been studied as a 
typical example to understand the common and generic 
feature. An important aspect of the transition is that 
it has both continuous and first-order transitions con- 
nected by the tricritical point (TCP). We have studied 
charge-order transitions in the extended Hubbard model. 
By introducing the onsite interaction U and the nearest- 
neighbor repulsion V, this model shows the transition 
to a charge-ordered phase at low temperatures with al- 
ternating charge density for the bipartite lattice, when V 



increases. In the classical limit of ignoring the itinerancy, 
the model has been solved first by the mean-field approx- 
imation and then on a square lattice by Monte-Carlo sim- 
ulation with the exchange Monte Carlo algorithm to get 
around the critical slowing down. 

Our results for critical exponents agree with those in 
the literature for the available studied quantities. On the 
lambda line, it reproduces the Ising exponents, r\ = 0.25, 
v = 1.0 and logarithmic divergence of the specific heat, 
which are consistent with Onsager's exact solution of the 
2D Ising model. Using Ehrenfest law, we have shown 
that this divergence of specific heat causes the logarith- 
mic divergence of the doublon susceptibility Xd(q)i the 
susceptibility of the nearest-neighbor charge correlation 
Xr and the charge susceptibility Xc on the lambda line. 
At the TCP, we obtain the order parameter exponents, 
r\ t = 0.11 ± 0.02, v t = 0.55 ± 0.01 while for the suscepti- 
bility of the order parameter exponent and the doublon 
susceptibility exponent defined in eq. (61) and eq. (62), 
we obtain 7t = 1.04 ± 0.06 and j 2t = 0.90 ± 0.06, respec- 
tively. Furthermore we clarify which fluctuations diverge 
at TCP both in the mean-field approximation and Monte 
Carlo calculations. For the mean-field approximation and 
the Monte Carlo calculations, we have shown that, the 
specific heat C, the susceptibility of the nearest-neighbor 
charge correlation xr and the charge susceptibility Xc di- 
verge with the same exponent as that of Xd(o)> namely 
72* at TCP. However Xc has a subtlety. Because of the 
particle-hole symmetry, Xc does not diverge at half filling, 
while, away from half filling, we confirm that Xc diverges 
at TCP. As far as the authors know, the divergence of 
Xc at TCP has not been recognized in the literature. 
The comparison of the mean-field results with the Monte 
Carlo results indicate that the both results show the di- 
vergences in the same quantities. Exact results by Monte 
Carlo calculations give only quantitative modifications of 
the critical exponents from the mean-field results. 

In the itinerant model, within the Hartree-Fock study, 
we obtain TCP at a finite temperature for t'/t = 0.2. 
At TCP, we obtain the critical exponent of the order 
parameter f3 t — 1/4. We show that the doublon sus- 
ceptibility Xd(o), and the susceptibility of the nearest- 
neighbor charge correlation xr diverge with the singu- 
larity <~ \g — g c \~ l2t at TCP. Here, 72* is equal to the clas- 
sical mean-field value 1/2. We also show that the charge 
susceptibility Xc diverges at TCP in the grand-canonical 
ensemble with the singularity ~ \g — g c |~ 72 where 72 = 1. 
The overall exponents are the same as those in the clas- 
sical model. 

At the TCP, various physical quantities diverge much 
more strongly than those at the lambda transition. 
In particular, diverging charge fluctuations at zero 
wavenumber revealed around TCP may mediate various 
instabilities toward ordering such as superconductivity, 
when the divergences are involved in the Fermi degen- 
eracy with Fermi surface instability. The consequence of 
such strong divergences clarified at the TCPs will be dis- 
cussed in a separate publication. 

We have also obtained the conductivity exponent, 
which is specific to the itinerant model. In three- 
dimensional systems, these exponents including the con- 



18 J. Phys. Soc. Jpn. 



Full Paper 



Online-Journal Subcommittee 



ductivity exponent may be correct because the upper 
critical dimension is three, while in two dimensions, the 
exponents may be modified as in the classical model, 
which we have shown in the Monte Carlo result. The con- 
ductivity exponent beyond the mean-field level in two- 
dimensional systems is left for future studies. 

On the conductivity, we obtain the exponent pt which 
is defined by \X - X c \ = \g - ff c |*(Alog \g - g c \ + B) 
in the canonical ensemble at TCP. The exponent p t has 
a relation to the order parameter exponent (3 t as p t = 
2p t = 1/2. We also show that p t changes into p' t = 2p t = 
1 at TCP in the grand-canonical ensemble when TCP in 
the canonical ensemble is in the phase separation region. 
Similarly, on the lambda line, the exponent defined by 
\X - X c \ = \g- g c \ p {A\og \g - g c \ + B) is given by p = 
2/3 = 1 on the lambda line, where p should be 2p t . 

The experimental data for (DI-DCNQI) 2 Ag 4 show q <~ 
1 for the points closest to the tricriticality, where q is de- 
fined from the conductivity a oc (g — g c ) q {A log \g — g c \ + 
B). There are two possibilities for explaining the exper- 
imental results. One is that the phase separation occurs 
and terminates at TCP in (DI-DCNQI) 2 Ag. In this case, 
the conductivity should have the same singularity as that 
of the Hartree-Fock analysis in the grand-canonical en- 
semble, namely, p' t = 1, which is consistent with the 
experimental value q ~ 1. Even when the long-ranged 
Coulomb repulsion prohibits the real phase separation, it 
is conceivable that the formation of microdomains results 
in effectively similar exponents in experiments. Another 
possibility is that the phase separation does not occur 
and terminates on the first-order transition line in (DI- 
DCNQI) 2 Ag. In this case, the conductivity should have 
the same singularity as that of the Hartree-Fock analysis 
in the canonical ensemble, namely p t = 1/2. The exper- 
imental results are not consistent with our Hartree-Fock 
value. Then the origin of this discrepancy must be as- 
cribed to the interpretation that the experimental data 
is not sufficiently close to the critical region of TCP. If 
one approaches the critical region of TCP, q should be 
close to pt — 1/2. To clarify the origin of singularity, it 
would be desired to perform more detailed experimental 
studies. 

We have discussed quantum effects through the sin- 
gularity of the conductivity, while the presence of the 
transfer does not alter the exponents for the order pa- 
rameter A itself from the classical value. This may be 
due to the fact that the transition temperature is suffi- 
ciently high so that the quantum proximity is not visi- 
ble for the charge-order transition itself. If the tricritical 
temperature could be lowered, its quantum effect would 
be observed. However, in the extended Hubbard model, 
the tricritical line terminates at a finite temperature and 
cannot be lowered to zero as one sees in the schematic 
phase diagram in Fig. 17. The quantum effect on TCP 
would be an intriguing issue to be studied in a different 
situation in the future. 

In the itinerant model, within the Hartree-Fock study, 
we have ignored the possibility of the magnetic order. 
However, antiferromagnetic and metal-insulator transi- 
tions of the extended Hubbard model with V = can 
also be studied in the framework by putting V = 0, where 



the interaction effect from U is only reflected in g in our 
treatments. Then the same argument for the criticali- 
ties apply in the region of g < 0, because transforming 
the charge order to the antiferromagnetic order with the 
transformation of g to —g leaves the self-consistent equa- 
tion unchanged. This means that the Hubbard model 
with only the onsite repulsion U with the transfers t 
and t' has the same phase diagram and criticalities sim- 
ply by replacing the charge order with the antiferromag- 
netic order. This is a one-to-one equivalence within the 
Hartree-Fock approximation. However, in the true phase 
diagram of the Hubbard model, it has a subtlety because 
the metal-insulator boundary (Mott transition) may ex- 
tend beyond the antiferromagnetic boundary contrary to 
the artifact of the Hartree-Fock results, while in case of 
the charge order, the insulating phase without the charge 
order should not exist and the Hartree-Fock phase dia- 
gram is correct in this aspect. 

We find that the tricritical line terminates when t' /t 
is increased and the phase boundary of metal-insulator 
transition separates from that of the charge-order tran- 
sition. This generates an intriguing first-order metal- 
insulator boundary with marginally quantum critical be- 
havior 14 within the mean-field level. This issue will be 
discussed in a separate paper. 
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Appendix A: Ehrenfest law 

We explain the general framework of Ehrenfest law 
in this appendix. 12 We consider the situation that the 
continuous transition line lies in a two-dimensional plane 
in a set of parameter space (see Fig. 32). We consider 
two physical properties X and Y, which arc conjugate 
to x and y, respectively. We define X and Y as the first 
derivatives of the free energy F with respect to x and y 
as 



x= dF y= 8F 

dx' dy 1 



Xx 



dX 

dx ' 



Xy 



dY 
dy ' 



(83) 



We assume that X and Y are continuous along the 
continuous transition line, but the derivatives of X and 
Y are not continuous. Then we define two limits of the 
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Fig. 32. Continuous transition line in two-dimensional plane. 



second derivatives as 
<9A, 



dx 

8X i 
dx 



\dis 



= lim 



dX 



x — >x~^ ,y — >y~^ 9x 



dx 

dx ' 



(84) 



where the subscripts dis and order express the deriva- 
tives in the disordered and ordered sides at the transition 
line, respectively. 

X is expanded by x and y near a critical point (x c , y c ) 

as 

X(x c + Sx, y c + Sy) = X(x c , y c ) 

+ ^r-\di S 5x + ^-\di S Sy + (8x 2 ,Sy 2 ) : disorder side, 
ax ay 

X(x c + Sx, y c + Sy) = X(x c , y c ) 



dX 

a^ 1 



r Sx + ^-\order5y + {Sx 2 , Sy 2 ) : order side. 
dy 

(85) 



Along the continuous line, X at the disorder and order 
phases have the same value at {x Cl y c ). We obtain the 
relation 

X(x + Sx,y + 5y)dis - X(x + Sx,y + Sy) order = 0. 

(86) 

From this relation, we obtain the slope of the continuous 
transition line. 

s ax I ax I 

Oy q x \dis q x \order 

Sx 



ax i _ ax i 

Q y \dis Q y \ order 
Xx\dis Xx\order 



DX 



DX I 



(87) 



I order 

Using the same discussion for Y, we obtain the relation 

as 

<■ §Y\ _ dY\ 

°V dx \ dls dx \order 

Sx 



dY | 
dy I 

dY | 

_ dx 



\order 



\dis 



dY | 



\order 



XY\dis ~ XY\order 



(88) 



Along the continuous "line" , second partial differen- 
tiations are well defined. Thus, the relations such that 

\dis , -Qy-\order — — \order are Satisfied. 



We assume that || is finite. In this condition, if xx 

diverges, ^ should diverge. Then, from the relation 

If = should diverge. Finally, \y also diverges 

with the same singularity of \x and vice versa. This is 
the Ehrenfest law. 

If we take X as the onsite interaction U and Y as the 
temperature T, xx and \y correspond to the doublon 
susceptibility xd(o) an d the specific heat C, respectively. 
From this, if the slope of the continuous transition line 
is finite, Xd(o) nas the same singularity as that of C. 

We summarize the relations with \x or xy and phys- 
ical properties below. 



XY = C (Y = T) 
Xx = Xd(o) (A = U) 
X x=XR (X = V) 
Xx =Xc (X = /i) 



(89) 
(90) 
(91) 
(92) 



Appendix B: Singularity of carrier density 

In this appendix, we obtain the singularity of the car- 
rier density X defined in eq. (77). We evaluate X near 
(7r/2,7r/2) and obtain logarithmic correction. 

Near (7r/2,7r/2), upper band E + is approximated as 

tt'SxSy + y/At 2 {8x + Sy) 2 + A 2 



E 4 



a 2 .) + v/4t 2 M 2 



A 2 , 



(93) 
(94) 



where wc define <5x[<5y] and ai[a2] as Sx = ir/2 — k x [Sy = 
ir/2 — k y ] and a.\ = Sx + Sy [012 = Sx — Sy], respectively. 

Using this band dispersion, the singular part of A is 
approximated as 

dct\doi2 



X 



-h 



+ cxp/3[t'(a 2 - a 2 ) + y / 4t 2 (ai) 2 + A 2 - fj] ' 

(95) 

Near the transition point, chemical potential can be ex- 
panded with respect to A as 

Mo + Mi A 2 , 



where /xo and fii are constants. 

To obtain a stronger singularity than that of A 2 , 
consider dX/dA 2 . 



(96) 



we 



ox 

dA 2 



= J da 1 da 2 g((3,E+) 



(97) 



\jAt 2 a\ + A 2 

where g{(3,E+) = -/3e^ E +-^/(l + e ^ E +-^f. Since 
g(fl,E + ) does not diverge at a finite temperature, its in- 
tegration should be O(l). Therefore, near the transition 
point (A <C 1), the strongest singularity of eq. (97) comes 
from / da l / y jAt 2 a\ + A 2 . Then we obtain the singular- 
ity of A as 

r A 



A oc 



/ 

Jo 



dai 



1 



yjAt 2 a{ + A 2 



log[«i - 
log[|A|/2i]/t, 



(A/2t) 2 ]/(2t) 



(98) 

(99) 
(100) 
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where A is a cutoff parameter. From this, we obtain the 
singularity of the carrier density X as 

A oc A 2 (AlogA + B), (101) 

where A and B are constants. 

Here, we note that, since the logarithmic correction 
comes from cq. (98), this logarithmic correction appears 
in any dimensions. 

Appendix C: Singularity of relaxation time 

In this appendix, we analyze the singularity of the co- 
herent part of the relaxation time r co h near TCP fol- 
lowing the conventional Fermi-liquid scheme. 1,18 In this 
scheme, although the singularity of r co h is caused by the 
antifcrromagnctic spin fluctuation, it is likely that the 
order-parameter fluctuation of the charge order causes a 
similar singularity. 

We now consider the singularity taking place from the 
disordered (non-charge-ordered) phase. The relaxation 
time T co h is proportional to the inverse of the imaginary 
part of self-energy ImS. In the Fermi-liquid scheme, ImE 
is approximated as 

\2 



imE~[c 2 + (^r) 2 ]]Tr Ti (k) 2 



(102) 



where Tf j(k) is the four-point vertex function. The four- 
point vertex function r^(k) has the relation with the 
susceptibility x(k) such that 

X(k) 



r u (k) 



(103) 



Xfc(k) 2 ' 

where Xfc(k) is a noninteracting susceptibility. Near the 
transition point, x(Q + q) can be approximated as 

1 



x(Q + q) 



(104) 



K 2 + q 2 ' 

Using this relation, we obtain the singularity of ImS as 

k d - 1 



Im£cx]rr n (k) s 



oc 



dk- 



(K 2 



oc K 



d-4 



(105) 



In a three-dimensional system, near TCP in the canonical 
ensemble, x(Q) diverges as 

x(Q)oc|T-T c |-T'. ( io6) 

From this, the singularity of K is given as 

K x\T -T c \^/ 2 , (107) 

where j t = 1- Therefore in a three-dimensional system, 
the singularity of r co h is obtained as 

rcoh oc T s coh = (ImE)- 1 oc \T - T c \ x ' 2 . (108) 

Here, if we employ the above analysis, T co h has no normal 
part r° oh . On the other hand, near TCP in the grand- 
canonical ensemble, j t — 1 changes into 7 = 2. There- 
fore, the singularity of r co h changes into 

r coh KT s coh = \T-T c \. (109) 



It should be noted that r co h has the same singularity 
as that of the carrier density X except logarithmic cor- 
rection in both the canonical and grand-canonical en- 
sembles. The singularity of X is logarithmically stronger 
than that of r co h- 

We, however, also note that the above analysis would 
oversimplify the real situation, where the momentum de- 
pendence of the relaxation time ignored in the above 
analysis could become important. In fact, if the pocket- 
type Fermi surface is expected in the ordered phase, 
strongly anisotropic arc-type structure may appear in the 
disordered phase. Then coherent quasiparticles on the arc 
along the Fermi surface cannot be scattered to other arc 
point by the scattering wave vector (n,ir). The diverg- 
ing fluctuation at (it, it) controls the quasiparticle relax- 
ation time, and the (tt, it) vector connects an arc point 
to the missing counterpart of the pocket, while such a 
missing part does not have a well-defined Fermi surface 
and such scattering process is not effective. If this cir- 
cumstance applies, r remains nonsingular and nonzero. 
Then we expect r to be expressed by r = r + t s , where 
r is the remaining regular part. We also expect a similar 
anisotropy in the self-energy in the ordered phase. 
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